/*
 * TotalEquilibriumNewtonProblem.h
 *
 *  Created on: 4 Sep 2011
 *      Author: Allan
 */

#ifndef TOTALEQUILIBRIUMNEWTONPROBLEM_H_
#define TOTALEQUILIBRIUMNEWTONPROBLEM_H_

// C++ includes
#include <vector>
using namespace std;

// Eigen includes
#include <Eigen/Core>
using namespace Eigen;

// GeoReact includes
#include "Core/Multiphase.h"
#include "Core/ReactionSystem.h"

#include "Numerics/NewtonSolver.h"

#include "Utils/Types.h"

// GeoReact forward declarations
class TotalEquilibriumSolver;

class TotalEquilibriumNewtonProblem : public NewtonProblem
{
public:
	TotalEquilibriumNewtonProblem();
	
	TotalEquilibriumNewtonProblem(const Multiphase& multiphase, const ReactionSystem& eReactions, const vector<string>& jSpecies, const vector<string>& eSpecies);
	
	virtual ~TotalEquilibriumNewtonProblem();
	
	void Initialiaze(const Multiphase& multiphase, const ReactionSystem& eReactions, const vector<string>& jSpecies, const vector<string>& eSpecies);
	
	void SetParameters(const VectorXd& uj, double T, double P, const VectorXd& n);
	
	virtual void Function(const VectorXd& nj, VectorXd& F);
	
	virtual void Jacobian(const VectorXd& nj, MatrixXd& J);
	
	void UpdateEquilibriumComposition(const VectorXd& nj);
	
	friend class TotalEquilibriumNewtonSolver;
	
private:
	/// The multiphase object
	Multiphase multiphase;
	
	/// The system of equilibrium-controlled reactions
	ReactionSystem eReactions;
	
	/// The number of primary (Nj) and equilibrium (Ne) species
	unsigned Nj, Ne;
	
	/// The indexes of the primary (jIndexes) and equilibrium (eIndexes) species
	vector<Index> jIndexes, eIndexes;
	
	/// The inverse of the stoichiometric matrix (uee)
	MatrixXd uee_inv;
	
	/// The canonical stoichiometric matrix (vej) = -inv(uee) * (uej)
	MatrixXd vej;
	
	/// The temperature (T) and pressure (P)
	double T, P;
	
	/// The mole vector (n) and activity vector (a) of all species
	VectorXd n, a;
	
	/// The equilibrium constants (Ke) and the reaction quotients (Qe) of the equilibrium-controlled reactions
	VectorXd Ke, Qe;
	
	/// The mole vector of the equilibrium species (ne)
	VectorXd ne;
	
	/// The total moles of the primary species
	VectorXd uj;
};

#endif /* TOTALEQUILIBRIUMNEWTONPROBLEM_H_ */
